#include <stdlib.h>
#include <cmath>
#include <vector>
#include <iostream>
using namespace std;

class Newton
{
 private:
  vector<double> x;
  vector<double> y;
  vector<vector<double>> table_dd;
  double t;
 
 public:
  double* answer;
 Newton(vector<double> _x, vector<double> _y)
    {
      x=_x;
      y=_y;
      int n=x.size();
      int i;
      int j;
      for (i=0; i<n; i++)
	{
	  vector<double> line;
	  for (j=0;j<n;j++)
	    {
	      line.push_back(0);
	    }
	  table_dd.push_back(line);
	}

      //The table of divided differences.
      for (i=0;i<n;i++)
	{
	  if (i==0)
	    {
	      for (j=0;j<n;j++)
		{
		  if (j==0)
		      table_dd[i][j]=y[i];
		  else 
		      table_dd[i][j]=0;
	        }
	    }
	  else
	    {
	      for (j=0;j<n;j++)
		{
		  if (j==0)
		      table_dd[i][j]=y[i];
		  else if (j>i)
		      table_dd[i][j]=0;
		  else
		      table_dd[i][j]=(table_dd[i][j-1]-table_dd[i-1][j-1])/(x[i]-x[i-j]);
	        }
	     }
	}
      for(i=0;i<n;i++)
	{
	  for(j=0;j<n;j++)
	    {
	      //cout<<"ssss"<<table_dd[i][j]<<endl;
	    }
	}
    
      double* p = new double[n];
      for (i=0;i<n;i++)
	{
	  if (i<n-1)
	    p[i]=0;
	}
      p[n-1]=table_dd[0][0];
      double* temp = new double[n];
       for (i=0;i<n;i++)
	{
	  if (i<n-2)
	    temp[i]=0;
	  else if (i==n-2)
	     temp[i]=1;
	  else
	    temp[i]=-x[0];
	}

       for (i=0;i<n-1;i++)
	 {
	   if(i>0)
	     {
	     double temp1[n];
	   for (j=0;j<n;j++)
	     {
	       temp1[j]=0;
	     }
	   double temp2[n];
	   for (j=0;j<n;j++)
	     {
	       temp1[j]=0;
	     }
	   //用temp进行系数计算
	   for (j=0;j<n-1;j++)
	     {
	       temp1[j]=temp[j+1];
	     }
	   
	   for (j=0;j<n;j++)
	     {
	       temp2[j]=temp[j]*(-x[i]);
	       // cout<<"temp1[j]"<<temp1[j]<<" temp2[j]"<<temp2[j]<<endl;
	       temp[j]=temp1[j]+temp2[j];
	     }
	   }
	   // for(j=0;j<n;j++)cout<<"temp "<<temp[j];
	   //cout<<endl;

	   for(j=0;j<n;j++)
	     {
	       //cout<<"ssss"<<p[i]<<" "<<temp[j]<<" "<<table_dd[i+1][i+1]<<endl;
	       p[j]=p[j]+temp[j]*table_dd[i+1][i+1];
	      
	     }
	 }

       answer=p;
    }

 
    double f(double _t)
  {
    t=_t;
    double dis=0;
    int n=x.size();
     for(int i=0;i<n;i++)
       {
         dis=dis+pow(t,n-1-i)*answer[i];
       }
    return dis;
  }
 
};


